Measurement of mitotic activity

ABSTRACT

A method of measurement of mitotic activity from histopathological specimen images initially identifies image pixels with luminances corresponding to mitotic figures and selects from them a reference pixel to provide a reference colour. Pixels similar to the reference colour are located; image regions are grown on located pixels by adding pixels satisfying thresholds of differences to background and image region luminances. Grown regions are thresholded in area, compactness, width/height ratio, luminance ratio to background and difference between areas grown with perturbed thresholds. Grown regions are counted as indicating mitotic figures by thresholding region number, area and luminance. An alternative method of measuring mitotic activity measures a profile of an image region and counts the image region as corresponding to a mitotic figure if its profile is above a threshold at an intensity associated with mitotic figures. A mitotic figure is also indicated if the profile does not meet the previous criterion but has three other values satisfying respective threshold criteria.

This invention relates to a method, an apparatus and a computer programfor measurement of mitotic activity, which indicates cell divisiontaking place in a tissue specimen: it is particularly relevant to makingmeasurements on potentially cancerous tissue such as breast cancertissue. The method is also relevant to other forms of cancer such ascolon and cervical cancer.

Breast cancer is a common form of female cancer, and also occurs to alesser extent in the male: once a lesion indicative of breast cancer hasbeen detected, tissue samples are taken and examined by ahistopathologist to establish a diagnosis, prognosis and treatment plan.However, pathological analysis of tissue samples is a time consuming andinaccurate process. It entails interpretation of images by human eye,which is highly subjective: it is characterised in particular byconsiderable subjectivity in observations of the same samples bydifferent observers and even by the same observer at different times.For example, two different observers assessing the same ten tissuesamples may give different opinions for three of the slides—30% error.The problem is exacerbated by heterogeneity, i.e. complexity of sometissue sample features.

Published International Application No. WO 02/47032 A1 relates tomeasurement of DNA in cells from cell images to indicate mitotic phase.It refers to the use of other image analysis parameters giving oneexample, magnitude of intensity variance, but does not give details ofhow this variance can be used to indicate mitotic activity.

There is a need to provide an objective form of measurement of mitoticactivity to inform a pathologist's diagnosis and patient treatment.

In one aspect, the present invention provides a method of measuringmitotic activity from histopathological specimen image data,characterised in that the method has the steps of:

-   -   a) identifying pixels in the image data having luminances        associated with mitotic figures;    -   b) selecting from among the identified pixels a reference pixel        which is sufficiently close in position and luminance to another        identified pixel to provide a reference colour;    -   c) locating pixels in the image data with luminances        sufficiently close to that of the reference colour to indicate        potentially mitotic figures;    -   d) incrementing image regions corresponding to potentially        mitotic figures from the located pixels by adding pixels        thereto, potential increments, to image regions being        implemented or rejected by according to whether or not their        luminances are sufficiently close to respective image region        luminances and sufficiently far from an image data background        luminance;    -   e) selecting grown image regions on the basis of thresholds for        image region area, compactness and width/height ratio; and    -   f) counting selected grown image regions as actually indicating        mitotic figures on the basis of a thresholds for number of such        regions.

The invention provides the advantage that it provides an objectivemeasurement of mitotic activity to inform a pathologist's diagnosis andpatient treatment.

The step of selecting grown image regions may also involve thresholdsfor ratio of image region luminance to background luminance and areadifference between areas derived by growing each image region withmultiple thresholds. The thresholds for image region area, compactness,width/height ratio, luminance and area difference may be: 355pixels<area<1700 pixels, 0.17<compactness<0.77, width/height ratio<2.7,luminance percentage<44%, area difference<23area/100.

The step of counting selected grown image regions may also involvethresholds for region area and luminance. Successive potentialincrements to image regions may be individual pixels each of which is animmediate row or column neighbour of an existing image region pixel.Step b) may implemented with, a reference pixel having a luminancediffering by less than 8% compared to another identified pixel distantfrom it by not more than two percent of a smaller of two imagedimensions.

Step a) may include white balancing and median filtering the image dataprior to identifying pixels having luminances corresponding to mitoticfigures. In step c) pixels may be cued for acceptance or rejection asregards indicating mitotic figures by:

-   -   a) thresholding colour image data to remove pixels lacking        intensities associated with mitotic figure imagery,    -   b) removal pixels not present in all colours, and    -   c) thresholding image region areas to remove those too small and        too large to be potential mitotic figures.

In step c) pixels may alternatively be cued for acceptance or rejectionas regards indicating mitotic figures by:

-   -   a) segmenting to identify pixels with intensities associated        with mitotic figure imagery,    -   b) thresholding image region areas to remove those too small and        too large to be potential mitotic figures,    -   c) cluster analysis to determine whether or not a pixel's image        region is in a sufficiently large cluster, and    -   d) necrotic and hairy edge filtering.

In another aspect, the present invention provides a method of measuringmitotic activity from histopathological specimen image data,characterised in that the method has the steps of:

-   -   a) measuring an intensity profile of an image region        corresponding to a potentially mitotic figure,    -   b) counting the image region as indicating a mitotic figure if        its profile has a value greater than a prearranged threshold at        a position in the profile having intensity associated with        mitotic figure imagery.

This aspect preferably includes counting the image region as indicatinga mitotic figure if its profile has a first value not greater than theprearranged threshold at a position in the profile having intensityassociated with mitotic figure imagery, a second value greater than aprearranged second threshold, a third value greater than a prearrangedthird threshold, and a minimum value less than a prearranged fourththreshold. The first value may be at one end of the profile, the firstand second values may adjoin one another in the profile and the thirdvalue may not adjoin the second value.

The image data may comprise a first Principal Component obtained byPrincipal Component Analysis (PCA) of coloured image data, and step a)may include preprocessing image data by:

-   -   a) decomposing the image data into overlapping sub-images,    -   b) applying PCA to the sub-images to derive a first Principal        Component image,    -   c) thresholding the first Principal Component image to produce a        binary image of blobs and background    -   d) rejecting blobs adjacent to or intersecting sub-image        boundaries,    -   e) filling holes in blobs,    -   f) rejecting blobs too small to correspond to potential mitotic        figures, and    -   g) reassembling the sub-images into a single image for image        region profile measurement as aforesaid in step a).

After step g) pixels may be cued for acceptance or rejection as regardsindicating mitotic figures by:

-   -   a) thresholding colour image data to remove pixels lacking        intensities associated with mitotic figure imagery,    -   b) removal pixels not present in all colours, and    -   c) thresholding image region areas to remove those too small and        too, large to be potential mitotic figures.

After step g) pixels may alternatively be cued for acceptance orrejection as regards indicating mitotic figures by:

-   -   a) segmenting to identify pixels with intensities associated        with mitotic figure imagery,    -   b) thresholding image region areas to remove those too small and        too large to be potential mitotic figures,    -   c) cluster analysis to determine whether or not a pixel's image        region is in a sufficiently large cluster, and    -   d) necrotic and hairy edge filtering.

In another aspect, the present invention provides computer apparatus formeasuring mitotic activity from histopathological specimen image data,characterised in that it is programmed to execute the steps of:

-   -   a) identifying pixels in the image data having luminances        associated with mitotic figures;    -   b) selecting from among the identified pixels a reference pixel        which is sufficiently close in position and luminance to another        identified pixel to provide a reference colour;!    -   c) locating pixels in the image data with luminances        sufficiently close to that of the reference colour to indicate        potentially mitotic figures;    -   d) incrementing image regions corresponding to potentially        mitotic figures from the located pixels by adding pixels        thereto, potential increments to image regions being implemented        or rejected by according to whether or not their luminances are        sufficiently close to respective image region luminances and        sufficiently far from an image data background luminance;    -   e) selecting grown image regions on the basis of thresholds for        image region area, compactness and width/height ratio; and    -   f) counting selected grown image regions as actually indicating        mitotic figures on the basis of a thresholds for number of such        regions.

Computer apparatus for measuring mitotic activity from histopathologicalspecimen image data, characterised in that it is programmed to executethe steps of:

-   -   a) measuring an intensity profile of an image region        corresponding to a potentially mitotic figure,    -   b) counting the image region as indicating a mitotic figure if        its profile has a value greater than a prearranged threshold at        a position in the profile having intensity associated with        mitotic figure imagery.

In yet another aspect, the present invention provides a computer programfor use in measuring mitotic activity from histopathological specimenimage data, characterised in that the computer program containsinstructions to control a computer to implement the steps of:

-   -   a) identifying pixels in the image data having luminances        associated with mitotic figures;    -   b) selecting from among the identified pixels a reference pixel        which is sufficiently close in position and luminance to another        identified pixel to provide a reference colour;    -   c) locating pixels in the image data with luminances        sufficiently close to that of the reference colour to indicate        potentially mitotic figures;    -   d) incrementing image regions corresponding to potentially        mitotic figures from the located pixels by adding pixels        thereto, potential increments to image regions being implemented        or rejected by according to whether or not their luminances are        sufficiently close to respective image region luminances and        sufficiently far from an image data background luminance;    -   e) selecting grown image regions on the basis of thresholds for        image region area, compactness and width/height ratio; and    -   f) counting selected grown image regions as actually indicating        mitotic figures on the basis of a thresholds for number of such        regions.

In an additional aspect, the present invention provides a computerprogram for use in measuring mitotic activity from histopathologicalspecimen image data, characterised in that its instructions provide forimplementing the steps of:

-   -   a) measuring an intensity profile of an image region        corresponding to a potentially mitotic figure, and    -   b) counting the image region as indicating a mitotic figure if        its profile has a value greater than a prearranged threshold at        a position in the profile having intensity associated with        mitotic figure imagery.

The computer apparatus and computer program aspects of the invention mayhave preferred features equivalent to corresponding method aspects ofthe invention.

In order that the invention might be more fully understood, embodimentsthereof will now be described, by way of example only, with reference tothe accompanying drawings, in which:

FIG. 1 is a block diagram of a procedure for measuring mitosis activityof the invention, and incorporating mitosis cueing, feature detectionand counting;

FIG. 2 shows in more detail mitosis cueing in the procedure of FIG. 1;

FIG. 3 is a block diagram of an alternative approach to mitosis cueingin the procedure of FIG. 1;

FIG. 4 illustrates the use of hidden Markov random field segmentation inthe FIG. 3 mitosis cueing procedure;

FIG. 5 shows in more detail mitosis feature detection in the procedureof FIG. 1; and

FIG. 6 is a block diagram of an alternative approach to mitosis featuredetection in the procedure of FIG. 1.

Referring to FIG. 1, a procedure 10 for the assessment of tissue samplesin the form of histopathological slides of potential carcinomas of thebreast is shown. This drawing illustrates processes for measuringmitotic activity to produce a parameter for use by a pathologist as thebasis for assessing patient diagnosis.

The procedure 10 employs a database 12, which maintains digitised imagedata obtained from histological slides as will be described later.Sections are taken (cut) from breast tissue samples (biopsies) andplaced on respective slides. Slides are stained using the staining agenthaematoxylin & eosin (H&E), which is a common stain for delineatingtissue and cellular structure. Tissue stained with H&E is used to assessmitotic activity.

Measurement of mitotic activity in a tissue specimen gives an indicationof the degree of cell division that is taking place. A histopathologicalslide is a snap shot representing a very short time interval in a celldivision process, so the chance of such a slide showing a particularphase of mitotic activity is very small: if such a phase is in factpresent in a slide, that is a good indicator of how fast a potentialtumour is growing.

In a prior art manual procedure for scoring mitotic activity, aclinician places a slide under a microscope and examines a region of it(referred to as a tile) at magnification of ×40 for indications ofmitotic activity. This manual procedure involves a pathologistsubjectively and separately estimating unusual colour, size, shape andboundary definition of cells in a tissue sample. The values obtained inthis way are combined by the pathologist to give a single measurementfor use in diagnosis. The process hereinafter described in this examplereplaces the prior art manual procedure with an objective procedure.

The invention uses image data from histological slides. In the presentexample, image data were obtained by a pathologist using Zeiss Axioskopmicroscope with a Jenoptiks Progres 3012 digital camera. Image data fromeach slide is a set of digital images obtained at a linear magnificationof 40 (i.e. 40×), each image being an electronic equivalent of a tile.

To select images, a pathologist scans the microscope over a slide, andat 40× magnification selects regions (tiles) of the slide which appearto be most promising in terms of analysing of mitotic activity. Each ofthese regions is then photographed using the microscope and digitalcamera referred to above, and this produces for each region a respectivedigitised image in three colours, i.e. red, green and blue (R, G & B).Three intensity values are obtained for each pixel in a pixel array toprovide an image as a combination of R, G and B image planes. This imagedata is stored temporarily at 12 for later use. Ten digitised images arerequired for measurement of mitotic activity at 14 which then providesinput to a diagnostic report at 28.

A number of alternative processes 16 to 24 will be described to measuremitotic activity in a given sample: these comprise two alternativemitotic cueing processes 16 and 18 and two alternative mitotic featuredetection processes 20 and 24. The measure of mitotic activity isconverted at 26 into a mitotic count for use by a pathologist.

Referring now to FIG. 2, the first alternative mitotic cueing process 16is shown in more detail. If it is selected for, use, it is carried outfor each of the ten digitised images mentioned above, but will now bedescribed for one such image referred to as the input image. It is usedto cue or identify dark image regions, which may correspond to mitoticcells.

At a stage 30, from the input image three histograms are formed showingthe occurrence frequency of pixel intensities, one histogramrepresenting R (red) intensities, one B (blue) and one G (green). Forexample, an image with 8 bits per colour per pixel would have ahistogram abscissa axis of 256 intensity values, 0 to 255, and ahistogram ordinate axis of number of pixels in the image having eachintensity value. Each histogram is a vector having 256 elements, and theith element (i=0 to 255) of each vector is the number of pixels havingintensity i in the R, G or B image plane.

The next stage is to threshold the R, G or B image planes at 32: toimplement this, firstly the total number N_(T) of pixels in an imageplane is counted (this will be the same value for all three imageplanes). For each image plane N_(T) is then divided by a respectiveempirical R, G or B parameter P_(R), P_(G) or P_(B) determined fromexperience of implementing the invention: parameter values P_(R)=100,P_(G)=100 and P_(B)=140 were derived manually and empirically from a setof 250 test images obtained using the Zeiss Axioskop microscope andJenoptiks Progres 3012 camera mentioned above. Images produced using adifferent microscope/camera combination might require differentparameters. This procedure gives three thresholds T_(R), T_(G) and T_(B)equal respectively to N_(T)/P_(R), N_(T)/P_(G) and N_(T)/P_(B).

The histograms and the thresholds T_(R), T_(G) and T_(B) are then usedfor each image plane to select low intensity pixels whose total numberdoes not exceed the threshold T_(R), T_(G) or T_(B). So for example animage having a total number of pixels N_(T) equal to 20,000 would have ared and green image planes with P_(R)=100, P_(G)=100 and T_(R) and T_(G)equal to 200. For the blue image plane T_(B) is 2×10⁴/140 or ˜142. In aneight-bit range of pixel intensities with values 0 to 255, the red imageplane histogram might have numbers of pixels 3, 20, 50, 7, 20, 80 and 65at pixel intensity values 0 to 6 respectively. The total number ofpixels having pixel intensity values 0 to 6 is 245, which exceeds thered image plane threshold T_(R) of 200; however the total number overpixel intensity values 0 to 5 is less than 200, and these are thereforeretained and pixels with intensity values 6 to 255 are rejected. Theprocedure retains a small part of the histogram, which corresponds tothe darker regions of the red image plane (mitotic cells tend to bedark). This procedure is repeated for the green and blue image planesusing their respective thresholds. The objective is to retain in eachimage plane a number of pixels which are likely to be in proportion tothe number of pixels in the image.

The next stage 34 is spatial filtering: here the red, green and blueretained pixels are compared and every pixel which is not retained inall three image planes after thresholding is rejected. Each pixelremaining after spatial filtering is then cued by assigning it a binary1 value and all other pixels in the image which have rejected are set tobinary 0: this creates a single combined binary image for output fromthe stage 34.

At 36, a technique known as “connected component labelling” (CCL) isapplied to the combined binary image from stage 34: this is a knownimage processing technique (sometimes referred to as ‘blob colouring’)published by Klette R., Zamperoniu P., ‘Handbook of Image ProcessingOperators’, John Wiley & Sons, 1996, and Rosenfeld A., Kak A. C.,‘Digital Picture Processing‘, vols. 1 & 2, Academic Press, New York,1982. CCL gives numerical labels to image regions which are “blobs” inthe binary image, blobs being groups of contiguous or connected pixelsof the same value 1 in a binary image containing 0s and 1s only: eachgroup or blob is assigned a number (label) different to those of othergroups to enable individual blobs to be distinguished from others. CCLalso provides blob areas in terms of number of pixels.

Blobs are then retained (pixels set to 1) or rejected (pixels set to 0)based on their dimensions: blobs are retained if they both contain from95 to 5000 pixels inclusive and have height and width of not more than2000 pixels. In this example the minimum area of 95 pixels isdeliberately set to a small value to avoid rejecting too many blobs thatmight be of interest for possible mitotic activity. The maximum area isset to a large value for the same reason. The output of stage 36 is abinary image containing a set of labelled blobs for analysis for mitoticactivity as will be described later. Stage 36 is useful for removingblobs which aren't likely to be relevant for later processing, but it isnot essential.

Referring now to FIG. 3, the second mitotic cueing process 18 is shownin more detail. As in the previous example, if the process 18 isselected for use, it is carried out for each of the ten digitised imagesmentioned above, but will now be described for one such image. At astage 40 a single iteration of a technique known as “Hidden MarkovRandom Field segmentation” is performed on the red component of anoriginal (RGB) input image of a tile. This segmentation is a known imageprocessing technique, see by Devijver P. A., ‘Image segmentation usingcausal Markov Random Field models’, Pattern Recognition, J.Kittler (ed)Lecture Notes in Computer Science 301, Springer-Verlag 1988; alsoDucksbury P. G., ‘Parallel model based segmentation using a 3^(rd) orderhidden Markov model’, 4^(th) IEE Int. Conf. Image Processing and itsApplications, Maastricht, 7-9 Apr. 1992. The input image is quantisedinto four levels (i.e. reduced from potentially 256 grey levels to justfour): the first of these levels corresponds to very dark areas of theinput image, second and third levels are progressively less dark and thefourth level corresponds to light image areas. The input image's greylevel histogram is computed and initially fitted with a set of fourGaussian distributions of number of pixels N(x) as a function of pixelvalue x, the distributions having a variety of means μ and standarddeviations a and each being of the form: $\begin{matrix}{{N(x)} = {\frac{1}{\sqrt{\left( {2\pi} \right)\sigma}}{\mathbb{e}}^{{{- {({x - \mu})}^{2}}/2}\sigma^{2}}}} & (1)\end{matrix}$

In this example, the Equation (1) distributions form the basis of imagesegmentation, which is defined as separation of objects from abackground in a digital image. The Gaussian distributions are numbered 1to 4 respectively. If an image pixel has a grey level that falls withinGaussian distribution number j (j=1, 2, 3 or 4), then its segmentedlabel is j. However where there is overlap of Gaussians as in shadedareas in FIG. 4 then the power of the Markov segmentation algorithmcomes into effect in making a probabilistic decision as to which is thebest segmented label to choose. The output from this stage 40 is a maskimage where each pixel has been allocated one of four numbers (1, 2, 3or 4). Only pixels having label 1 (corresponding to dark image areas)are used in subsequent processing stages. Even though only label 1pixels are retained, for good image segmentation it is better to use anadequate number of levels or labels (four in this case) and to assignpixels correctly: this improves the results for label 1. Using too fewsegmentation levels, segmented regions may be too large and containinappropriate pixels; with too many levels, segmented regions may be toosmall and fragmented.

Segmentation into four is found to be a good compromise, despite usingonly lowest numbered pixels in later processing. Stage 40 is a usefulsegmentation technique but others can be used instead.

At 42 a K-means clustering process is performed: this is a known processwhich is described by J. A. Hartigan and M. A. Wong in a paper entitled‘A K-means clustering algorithm’, Algorithm AS 136, Applied StatisticsJournal, 1979. K-means is an iterative statistical technique forcomputing an optimal set of clusters (groups of data items sharing somecommon property) from a dataset. This process uses raw image pixels fromthe red component of the original (RGB) input image of a tile input at40: it selects those red pixels that are located in the same imageposition as pixels labelled 1 in the mask image generated at 40.Designating the maximum and minimum of the values of the label 1 rawimage pixels as “mae” and “min” respectively, three clusters are usedand initially cluster centres for the K-means technique are set relativeto the dynamic range (max-min) of these pixels as follows: ClusterInitial cluster centre setting 1 min + 0.1 × (max − min) 10% of dynamicrange above min 2 min + 0.5 × (max − min) middle of dynamic range 3 max− 0.1 × (max − min) 10% of dynamic range below max

From the result of the K-means technique-two masks are created: one maskmarks areas containing cluster 1 OR cluster 2 (the joint cluster mask)whilst the second mask marks just cluster 1 in isolation (the cluster 1mask). Only clusters 1 and 2 are required for later processing; thethird cluster is used simply to separate data on a more reasonablebasis.

At 44 connected component labelling (CCL, as previously described) isapplied to the mask created at 42 which marks the joint cluster (cluster1 OR cluster 2). CCL also gives the areas of labelled blobs in terms ofnumbers of pixels. Blobs that are outside an allowed area range arerejected, and those within this range are retained as indicated in thetable below and renumbered sequentially from label 1 upwards. The outputfrom stage 44 is a set of labelled blobs. EITHER blob size < 70 Rejectblob - set all its pixels to 0 and pixels OR size > 262144 exclude fromset of labelled blobs pixels 70 pixels ≦ blob Accept for furtherprocessing and include size ≦ 262144 pixels in set of labelled blobs

The maximum and minimum acceptable blob areas (262144 and 70 pixels)were chosen to be sufficiently widely separated to avoid significantloss of potentially relevant blobs—this dynamic range might be reducibleif more knowledge of relevant blob sizes becomes available. Stage 44 isuseful to remove unwanted blobs but could be dispensed with.

At 46 each labelled blob generated at 44 is considered and accepted orrejected based on a contextual analysis of its cluster 1/(cluster1+cluster 2) score; i.e. from the results of 42, for each blob find howMany of its pixels are in cluster 1 (NC1) and how many in cluster 2(NC2): then calculate the ratio NC1/(NC1+NC2), and accept the blob forfurther processing if this ratio is 0.6 or greater. Reject the blob ifthe ratio is less than 0.6. The result is a reduced set of blobs.Contextual measure ≧ cluster context Accept for further processingthreshold (60%) into set of blobs Contextual measure < cluster contextReject blob threshold (60%)

At 48 take the reduced set of blobs generated at 46 and perform arejection using a ‘necrotic filter’: such a filter computes a standardmetric referred to as the Euclidean metric M_(Eu) of normalisedquantised boundary phase (as described later). The Euclidean metricM_(Eu) is of the form: $\begin{matrix}{M_{Eu} = \sqrt{\sum\limits_{i = 1}^{n}\left( {x_{i} - y_{i}} \right)^{2}}} & (2)\end{matrix}$where x_(i) and y_(l) (i=1 to n) are elements of vectors x and yrespectively representing the two quantities to be compared. The‘necrotic filter’ process is as follows: firstly, a Sobel edge filter isapplied to the labelled image obtained from step 44 and the rawimage—i.e. the red component of the original (RGB) input image input tostage 40 tile. The labelled image is used to obtain the boundary of theblob that may correspond to a mitotic figure and the raw image is usedto obtain the phase angle of pixels in the raw image. Sobel is astandard image processing technique published in Klette R., & ZamperoniP., ‘Handbook of image processing operators’, John Wiley & Sons, 1995. ASobel filter consists of two 3×3 arrays of numbers S_(P) and S_(Q), eachof which is convolved with successive 3×3 arrays of pixels in an image.Here $\begin{matrix}{S_{P} = {{\begin{bmatrix}1 & 2 & 1 \\0 & 0 & 0 \\{- 1} & {- 2} & {- 1}\end{bmatrix}\quad{and}\quad S_{Q}} = \begin{bmatrix}1 & 0 & {- 1} \\2 & 0 & {- 2} \\1 & 0 & {- 1}\end{bmatrix}}} & (3)\end{matrix}$

A first 3×3 array of pixels is selected in the top left hand corner ofthe labeled image:

designating as C_(ij) a general labelled pixel in row i and column j,the top left hand corner of the image consists of pixels C₁₁ to C₁₃, C₂₁to C₂₃ and C₃₁ to C₃₃. C_(ij) is then multiplied by the respective digitof S_(P) located in the S_(P) array as C_(ij) is in the 3×3 cyan pixelarray: i.e. C₁₁ to C₁₃ are multiplied by 1, 2 and 1 respectively, C₂₁ toC₂₃ by zeroes and C₃₁ to C₃₃ by —1, −2 and −1 respectively. The productsso formed are added algebraically and provide a value p. The value of pwill be relatively low for pixel values changing slowly between thefirst and third rows either side of the row of C₂₂, and relatively highfor pixel values changing rapidly between those rows: in consequence pprovides an indication of edge sharpness across rows. This procedure isrepeated using the same pixel array but with S_(Q) replacing S_(P), anda value q is obtained: q is relatively low for pixel values changingslowly between the first and third columns either side of the column ofC₂₂, and relatively high for pixel values changing rapidly between thosecolumns, and q therefore provides an indication of edge sharpness acrosscolumns. The square root of the sum of the squares of p and q are thencomputed i.e. √{square root over (p²+q²)}, which is defined as an “edgemagnitude” and becomes T₂₂ (replacing pixel C₂₂ at the centre of the 3×3array) in the transformed image. Tan⁻¹p/q is also obtained at each pixeland is defined as a “phase angle”.

A general pixel T_(ij) (row i, column j) in the transformed image isderived from C_(i−1,j−1) to C_(i−1,j+1), C_(i,j−1) to C_(i,j+1) andC_(i+1,j−1) to C_(i+1,j+1) of the labeled image. Because the central rowand column of the Sobel filters in Equation (3) respectively are zeros,and other coefficients are 1s and 2s, p and q for T_(ij) can becalculated as follows:p={C _(i−1,j−1)+2C _(i−1,j) +C _(i−1,j+1) }−{C _(i+1,j−1)+2C _(i+1,j) +C_(i+1,j+1)}   (4)q={C _(i−1,j−1)+2C _(i,j−1) +C _(i+1,j−1) }−{C _(i−1,j+1)+2C _(i,j+1) +C_(i+1,j+1)}   (5)

Beginning with i=j=2, p and q are calculated for successive 3×3 pixelarrays by incrementing j by 1 and evaluating Equations (2) and (3) foreach such array until the end of a row is reached; j is then incrementedby 1 and the procedure is repeated for a second row and so on until thewhole image has been transformed. The Sobel filter cannot calculatevalues for pixels at image edges having no adjacent pixels on one orother of its sides: i.e. in a pixel array having N rows and M columns,edge pixels are the top and bottom rows and the first and last columns,or in the transformed image pixels T₁₁ to T_(1M), T_(N1) to T_(NM), T₁₁to T_(1M) and T_(1M) to T_(NM). By convention in Sobel filtering theseedge pixels are set to zero. The output of the Sobel filter comprisestwo transformed images, one (the edge filtered image) contains theboundaries of the labeled blobs produced at 44 whilst the other containsthe “phase angle” of the raw input image.

A pixel of a labelled blob is a boundary pixel if the like-located pixelin the Sobel edge filtered image is non-zero: for each boundary pixelthe phase angle is extracted from the like-located pixel in the Sobelphase angle image. This phase angle information is then quantised toreduce it to four orientation ranges (0-44, 45-89, 90-134, 135-179degrees) and the number of boundary pixels in each orientation isnormalised by dividing it by the number of pixels in the perimeter ofthe blob. This results in a respective 4-element vector of normalisedorientations or quantised phase for each blob: The Euclidean measure ofeach of these vectors is computed using Equation (2) and compared withthat of a perfect circle, in which the four vector elements are of equalvalue. This searches for blobs relatively far from circularity. TheEuclidean measures of the labelled blobs are compared with a EuclideanThreshold of 0.4354 and are rejected if they are greater than it. TheEuclidean threshold was derived from a K-means analysis of a testdataset. Output a set of remaining blobs. Euclidean measure > 0.4354Reject blob Euclidean measure ≦ 0.4354 Accept blob for furtherprocessing as part of a set of remaining blobs

At 50 each of the blobs remaining after step 48 is examined and a ‘hairyedge filter’ operation is performed. A ‘hairy edge filter’ measures theamount of edge structure in an area around a blob, this being a roughapproximation to the ‘hairy fibres’ sometimes seen around mitoticfigures. This is computed for each blob as follows:

-   -   a) Two morphological dilations are applied using filters 5×5 and        13×13 pixels in size as shown below. $\begin{matrix}        {{{5 \times 5}\quad{filter}} = \begin{matrix}        0 & 0 & 1 & 0 & 0 \\        0 & 1 & 1 & 1 & 0 \\        1 & 1 & 1 & 1 & 1 \\        1 & 1 & 1 & 1 & 0 \\        0 & 0 & 1 & 0 & 0        \end{matrix}} & (6) \\        {{{13 \times 13}\quad{filter}} = \begin{matrix}        0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\        0 & 0 & 0 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 0 & 0 & 0 \\        0 & 0 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 0 & 0 \\        0 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 0 \\        0 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 0 \\        0 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 0 \\        1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 \\        0 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 0 \\        0 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 0 \\        0 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 0 \\        0 & 0 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 0 & 0 \\        0 & 0 & 0 & 1 & 1 & 1 & 1 & 1 & 1 & 1 & 0 & 0 & 0 \\        0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0        \end{matrix}} & (7)        \end{matrix}$

Morphological dilation is an expansion operation: for an original binaryimage (i.e. having pixel values 1 and 0 only), the expansion operationcomprises locating each pixel with value 1 and setting pixels in itsvicinity also to 1. In the two above arrays, a central 1 indicates apixel found to be 1 in the image, other 1s indicate the relativepositions of nearby pixels set to 1 to implement morphological dilation,and 0s represent pixels left unchanged.

Morphology is an image processing technique based on shape and geometry.It is a standard image processing procedure published in Umbaugh S. C.,‘Colour vision and image processing’, Prentice Hall, 1998. Morphologyapplies a filter of some size and shape to an image. In the simplestsense dilation (dilates (or expands) an object) at each pixel positionthe output of the dilation is the logical OR of the inputs. The filtersused contain approximations to circles as shown above. The applicationof the two dilation operations resuits in two dilated results. Each blobis dilated by two different amounts, as described with reference toEquations (6) and (7): the blob resulting from the 5×5 filter is thensubtracted from that resulting from the 13×13 filter, which results in aborder around the blob. This is repeated for each of the blobs remainingafter step 48.

b) A Sobel filter (as previously described) is applied to obtain thegradient of the raw image—i.e. the red component of the original (RGB)input image of a tile input at 40. A summation of gradient values isformed within the area of the border (mask) determined at a) around eachblob and is used as a filter measure for comparison with a HairyThreshold for acceptance and rejection of blobs as shown in the tableimmediately below. Filter measure ≧ Hairy Threshold Accept blob as validfor processing (25000) for mitotic feature detection Filter measure <Hairy Threshold Reject blob (25000)

The output of step b) is a set of cued blobs considered valid for use inmitosis feature detection as will be described later. Stages 46, 48 and50 are desirable to reduce unwanted blobs, but are not essential if theconsequent processing burden in mitotic feature detection can betolerated.

Referring now to FIG. 5, there is shown a flow diagram of the mitoticfeature detection process 20, which is carried out for each of the tendigitised input images referred to above and will be described for oneimage. At a first stage 60, an input RGB image is preferably whitebalanced by remapping its most luminous pixel to white. For each pixel i(i=1 to total number of pixels in image) in this image, pixel luminanceL_(i) is computed from its red, green and blue Intensities R_(i), G_(i)and B_(i) using the following equation:L _(i)=0.299×R _(i)+0.587×G _(i)+0.114×B _(i)   (8)

Then the pixel with the maximum of the luminance values of all pixels inthe input (RGB) image is located and used to record the correspondingvalues of R, G and B at that maximum luminance pixel position anddenoted LumMaxR, LumMaxG and LumMaxB. The ratios for each of the threeimage planes are then computed asRatioR=(255/LumMaxR)×1.05   (9)RatioG=(255/LumMaxG)×1.05   (10)RatioB=(255/LumMaxB)×1.05   (11)

The original RGB pixel values are now multiplied by these ratios toproduce a white balanced image with three image planes with thefollowing values for each pixel i:BalancedR _(i) =R _(i)×RatioR   (12)BalancedG _(i) =G _(i)×RatioG   (13)BalancedB _(i) =B _(i)×RatioB   (14)

The final stage is to clip the new white balanced image so that no pixelvalues lie outside the eight bit range (0 to 255). If any pixel value isless than 0, it is set to zero and if any pixel value is greater than255, it is set to,255. Production of a white balanced image is notessential but desirable to reduce variation between images.

At 62 the clipped white balanced image from step 60 is filtered with a3×3 median filter to remove spatial noise (desirable but not essential).The filter is applied independently to each of the balanced red(BalancedR), green (BalancedG) and blue (BalancedB) image planescomputed at 60: the median filter operation selects each pixel in theseimage planes in succession (other than edge pixels) and takes a 3×3array of those pixels centred on the selected pixel. The 3×3 array ofpixels is then sorted into ascending order of pixel value using what isreferred to as “quicksort”. Quicksort is a known technique published byKlette R., Zamperoniu P., ‘Handbook of Image Processing Operators’, JohnWiley & Sons, 1996, and will not be described. It is not essential butconvenient. The median pixel value (fifth of nine) is then taken as thefilter output to replace the value of the selected pixel. This isrepeated across the clipped white balanced image. Pixels in edge rowsand columns do not have the requisite 3×3 array, and for these theclipped white balanced image pixel values are retained in the medianfiltered image.

At 64 an ‘autopick colour’ process is applied which picks or locates apixel having the lowest luminance (darkest) in the median filtered image(excluding outlying pixels relatively remote from pixels of similarluminance): this means that the chosen pixel has at least one relativelynearby pixel with luminance similar to its own. Dark pixels are selectedbecause mitotic figures tend to have relatively low luminance usingconventional histological slide preparation techniques. The computationis as follows, for each pixel position in the median filtered image theluminance L_(i) is computed as follows:L _(i)=0.299×R _(i)+0.587×G _(i)+0.114×B _(i)   (15)

Computing L_(i) for each pixel position in the median filtered imageusing Equation (15) provides a luminance image: in the luminance image afirst pixel value and its location is stored as a current darkest pixel.Successive pixels in that image are compared with the first pixel: ifany comparison pixel has L_(i) darker (lower value of luminance) thanthe current darkest then its pixel value and its location are stored ina list of darkest pixels. After the list has reached ten entries, thecurrent least dark pixel in it is removed and replaced by a later pixelon each occasion a later pixel is darker than the least dark pixel: thisprocess continues until all pixels have been compared with and whereappropriate added to the list. After processing the entire luminanceimage, the list of darkest pixels is sorted into descending order ofdarkness using Quicksort (as mentioned above) so that the darkest isfirst. If the procedure results in, less than ten darkest pixels, thestored comparison luminance value originally obtained from the firstpixel is increased and the procedure repeated.

The next step is to determine whether or not the darkest pixel satisfiesthe condition that it is located relatively near another of the tendarkest pixels: i.e. the condition is that these two pixels be separatedby a distance of not more than twenty pixels in any direction (along arow, a column diagonally or between a diagonal and a row or column).This condition applies to an image of dimensions 1476 pixels by 1160pixels, so the maximum separation is 2% of the smaller image dimension.If this condition is satisfied, a darkest pixel has been located havinga comparatively near neighbour of similar luminance and therefore notconsidered to be an outlying pixel: the luminance of this darkest pixelis denoted by L_(picked colour). If the condition is not satisfied thenthe procedure is iterated by discarding the current darkest pixel fromthe list and taking the remaining darkest pixel; iteration continuesuntil the condition is met. The list size of ten was determined byanalysis of the ten images selected as previously described, but it isnot critical and a different number can be chosen. The pixel from whichL_(picked colour) is taken should therefore be selected from a smallgroup (twenty or less) of the darkest pixels.

At 66 a ‘colour proximity highlighting’ is carried out which locatesimage pixels in the median filtered image that have luminances differingfrom L_(picked colour) by less than 20 (i.e. less than 8% for aneight-bit intensity range from 0 to 255). This is carried out bycreating a mask image as follows: for each pixel in the median filteredimage, if a luminance inequality condition (16) below is true then thepixel is accepted and represented by a value 255 in the equivalentposition in the mask image.|(0.299×R _(i)+0.587×G _(i)+0.114×B _(i))−(L _(picked colour))|<20  (16)

If the inequality condition (16) is not true, then the current pixel isrejected and represented by 0 in the mask image. Results from whicheverof the cueing processes 16 and 18 has been used are introduced in thisstep 66: i.e. a pixel for which the inequality condition (16) is true isaccepted or rejected according to whether there is a 1 or 0 valuerespectively for a corresponding pixel located in the same position inthe mitosis cued image resulting from process 16 or 18. It is notessential to use a process 16 or 18 but it is useful to reduceprocessing burden.

At 68 accepted pixels are “grown” so that they come to correspond to awhole cell instead of just part of a cell. Luminance proximity is usedto check if growing should continue: i.e. if the mask image computed at66 indicates that there is too large a luminance difference between aselected pixel and a test pixel as compared to a luminance differencethreshold denoted by L_(T) (L_(T)=75 in this example), growth with theselected pixel does not continue.

The process 68 of growing pixels is as follows: firstly, an image storelabelled ‘grow’ is created that indicates whether pixel positions are‘grown’ or ‘not grown’ and each entry or pixel in ‘grow’ is initiallyset to ‘false’ (false=0) indicating that no pixels have yet been‘grown’: ‘true’ (true=1) would indicate ‘grown’. Secondly, a backgroundcolour for image ‘grow’ is computed from the median filtered image byaveraging all pixels other than those that are white (i.e. having R, Gand B all equal to the maximum value of 255): the value of thisbackground colour is recorded. The growing process now proceeds inaccordance with the computer program steps below: in these steps aconvention is used that an inset of a line to the right indicates aniterative loop including the line and those following it of equal andgreater inset, the loop terminating when a line of lesser inset isreached.

For each pixel in the median filtered image:

-   -   If the mask image pixel is accepted and a corresponding pixel in        the same location in ‘grow’ is 0 then:        -   Set the pixel in ‘grow’ to true (1)        -   Enter the position of the pixel in ‘grow’ in a list entitled            ‘Action List’        -   As long as the Action List is not empty then:            -   Remove the most recently added pixel from the Action                List but retain its image position in memory and                designate it as the ‘removed pixel’,            -   Select four pixels in ‘grow’ that are nearest neighbours                of the removed pixel, i.e. pixels immediately adjacent                to the removed pixel in the same column or row only, not                diagonal neighbours.                -   Select one of the nearest neighbour pixels not yet                    compared with the mask image pixel,                -   If (and only if) the three following criteria are                    met, i.e. (a) the mask image pixel and the selected                    nearest neighbour pixel differ in luminance by less                    than L_(T), (b) the mask image pixel and the                    background colour differ in luminance by more than                    L_(T), and (c) the mask image pixel and a white                    pixel differ in luminance by more than L_(T), then:                -    Change from false (0) to true (1) the pixel in                    ‘grow’ located in the same position as the selected                    nearest neighbour pixel, and add the position of the                    selected nearest neighbour pixel to the Action List                -   If any one or more of the criteria (a), (b) and (c)                    are not met, leave the relevant pixel in ‘grow’                    unchanged and do not add its position to the Action                    List                -   Repeat for next nearest neighbour pixel not yet                    compared with the mask image pixel,            -   Repeat for next entry in the Action list        -   Repeat for all other mask image pixels.

The above computer program steps provide a mechanism for continuing togrow a cell from an original single pixel in the mask image byreassessing further pixels for growth. Nearest neighbour pixels of a‘seed’ mask image pixel to be grown are assessed: each of the nearestneighbours which becomes added (changed to 1) in growth is also added tothe Action List for its uncompared nearest neighbours to be assessed.Growth therefore proceeds until all pixels adjoining but not part of agrown cell have been assessed and have failed one or more of the threeluminance criteria (a), (b) and (c). Growth then terminates for thatcell and restarts for another cell based on a new ‘seed’ mask imagepixel.

The result of 68 is a new image ‘grow’ which now contains a set of blobs(image regions of contiguous pixels of value 1) which are candidates forindicating positions of real cells that are likely to be of interest formitosis. The blobs are processed at 70 by connected component labellingas previously described: this derives a set of measurements for eachblob as follows:

-   -   area A (the number of pixels in a blob),    -   perimeter P (the number of pixels on a blob's boundary),    -   compactness (4λA/P²),    -   width (maximum number of pixels in a row across a blob),    -   height (maximum number of pixels in a column down a blob),    -   ratio width/height,    -   luminance percentage: using median filtered image pixels located        in the same positions as the pixels of the blob, this measure is        computed by multiplying by 100 the result of dividing the        luminance of the darkest median filtered image pixels by the        luminance of the background colour,    -   perturbed difference (difference between grown blob sizes        obtained using thresholds corresponding to L_(T) perturbed by a        prearranged increment and decrement respectively), and    -   Hue difference (absolute value of difference between an average        Hue and a background Hue from the input image).

The perturbed difference is computed as follows:

-   -   the threshold L_(T) is adjusted by adding a perturbation factor        P_(F) (P_(F)=4) in this example,    -   the growing process 68 is applied resulting in a new larger        blob,    -   the threshold L_(T) is adjusted by subtracting the perturbation        factor P_(F), and    -   the growing process 68 is applied resulting in a new smaller        blob,

A logical EXOR function is then computed between the new larger andsmaller blobs: i.e. each pixel in the smaller blob is EXORed with arespective pixel in the same position in the larger blob. Outer pixelsof the larger blob for which there are no like-located pixels in thesmall blob are treated as being EXORed with a different pixel value. TheEXOR function yields a 1 for a pair of pixels of different value and a 0for a pair of pixels of the same value. Its results provide an EXORimage with each EXOR value located as a pixel in the same position asthe blob pixels giving rise to it. The number of pixels equal to 1 inthe EXOR image is then counted and this number is the perturbeddifference.

The Hue difference is obtained for each blob as follows: the averagecolour of the median filtered image pixels located in the same positionsas the pixels of the blob is computed. This average colour and thebackground pixel colour obtained earlier are then converted fromred/green/blue (RGB) to a different image space hue/saturation/value(HSV). The RGB to HSV transformation is described by K. Jack in ‘VideoDemystified’, 2^(nd) ed., HighText Publications, San Diego, 1996. Inthis example the V and S components are not required. H is calculatedfor the average colour of each blob and the background pixel colour asfollows:Let M=maximum of (R,G,B)   (17)Let m=minimum of (R,G,B)   (18)Thennewr=(M−R)/(M−m)   (19)newg=(M−G)/(M−m)   (20)newb=(M−B)/(M−m)   (21)

Hue (H) then given by:If R equals M then H=60(newb−newg)   (22)If G equals M then H=60(2+newr−newb)   (23)If B equals M then H=60(4+newg−newr)   (24)If H greater than or equal 360 then H=H−360   (25)If H less than 0 then H=H+360   (26)

The difference between the H values of the average colour of the medianfiltered image pixels in the blob and the background pixel is thencalculated for each blob and becomes designated as the Hue differencefor that blob.

If a blob's parameters satisfy all the conditions in the tableimmediately below then the blob is accepted for further processing,otherwise it is rejected (deleted) by setting all its pixels to 0.Parameters for Accepted Blobs 355 pixels < Area < 1700 pixels 0.17 <Compactness < 0.77 Ratio of width/height < 2.7 Luminance percentage < 44Perturbed difference < Area × 23/100 Hue difference > 0

In the present example, the Hue difference will virtually always be true(due to the zero threshold in the table immediately above). However insome circumstances it may be desirable to have a non-zero threshold.Luminance percentage, perturbed difference and Hue difference are notessential, and can be omitted from the thresholds in the table abovegoverning further processing or otherwise.

At 72 a two quicksorts (as previously defined) are applied to the blobsto sort them into two lists, one of blobs in ascending order of blobarea and the other of blobs in ascending order of blob luminance in themedian filtered image. The blobs in the two lists are now referred to inaccordance with groupings that they mark the end of: i.e. a blob isreferred to as the “darkest X % blob” to indicate that it together withblobs (if any) of lower luminance than it are X % of the total number ofblobs. Similarly, a blob is referred to as the “largest Y % blob” toindicate that it together with blobs (if any) of greater area than itare Y % of the total number of blobs.

A blob of median area (the “median blob”, central in the area list) isnow identified. If there is an even number of blobs in the area list,the average area of the two central blobs is taken as the median area.Also at 72, unwanted blobs are eliminated to leave those assessed ascorresponding to mitotic figures as follows: if more than A blobs arepresent and the largest blob is more than B percent of the area of themedian blob, we retain each of the largest C blobs which has a luminancenot greater than that of the darkest D percent blob. Otherwise, if thelargest blob is not more than B percent of the area of the median blob,retain each of the darkest E blobs which has an area less than or equalthe largest F percent blob. In this example, values for A, B, C, D, Eand F are A=2 blobs, B=200 percent, C=3 blobs, D=30 percent, E=2 blobs,F=80 percent. The process 72 is computed as follows:

If the number of blobs is less than or equal to A then accept each blobas being a mitotic figure, the mitotic figure count is A and processingthe current image terminates.

Otherwise, if the number of blobs is greater than A:

-   -   obtain area of largest blob,    -   obtain area of median blob,    -   if the largest blob area is more than B percent of the median        blob area, then:        -   obtain the luminance of the darkest D percent blob—if there            is no blob at a        -   D percent position take the luminance of a blob which is            both darker and nearest to a notional D percent blob            position in the luminance list,        -   for each of the largest C blobs,            -   if blob luminance is less than or equal to the darkest D                percent blob,                -   accept the blob as being a detected mitotic figure,                -   increase a count of mitotic figures by 1, and                -   repeat for remainder of largest C blobs until none                    remain unassessed and output count of mitotic                    figures.

Otherwise, if largest blob area is not more than B percent of medianblob area,

-   -   obtain the area of the largest Fpercent blob—if there is no blob        at an F percent position take the area of a blob which is        nearest (in area) to a notional Fpercent blob position in the        area list,    -   For each of the darkest E blobs        -   if blob area is less than or equal to area of the largest F            percent blob,            -   accept the blob as being a detected mitotic figure, and            -   increase a count of mitotic figures (initially 0) by 1,            -   repeat for remainder of the largest E blobs until none                remain unassessed and output count of mitotic figures.

The criterion number of blobs is not greater than A may be the only oneused if desired, the mitotic count being taken as zero if this is notsatisfied. The “otherwise” criterion, i.e. number of blobs greater thanA and subsequent criteria, provide a further option.

As previously mentioned, the process 20 is carried out for a total often images or tiles: this repetition is to increase the likelihood ofobserving mitotic activity. The mitotic counts for the ten images arethen added together to provide a sum which is converted to a mitoticactivity grading as will be described later.

Referring now to FIG. 6, there is shown a flow diagram of an alternativemitotic feature detection process 24 which is carried out for each ofthe ten digitised images mentioned above. Although it is possible to usethe whole image, at 120 the digitised image (hereinafter the “inputimage”) is for convenience separated into overlayping windows of size128×128 pixels. The windows overlap with 64 pixels in both horizontaland vertical directions. Thus, each window overlaps half of itspreceding window above and to the left (if available). In each windowPrincipal Component Analysis (PCA, Karhunen-Loeve Transform) is applied.PCA is a standard mathematical technique described by Jollie I. T.,‘Principal Component Analysis’, Springer series in statistics, SpringerVerlag, 1986. It is also described by Jackson J. E., ‘A User Guide toPrincipal Components’ pp 1-25, John Wiley & Sons, 1991. PCA is atechnique for transforming a set of (possibly correlated) variables intoa smaller number of uncorrelated variables called principal components.The first principal component accounts for as much of the variability inthe set of variables as possible as compared to the other components:because of this it can be superior to taking a red, green or blue imageplane or an average thereof, which are also options at this stage. PCAinvolves obtaining a covariance matrix of the set of variables andsolving for its eigenvalues and eigenvectors. The covariance matrix iscalculated using the formula $\begin{matrix}{C_{i,j} = \frac{\sum\limits_{k = 1}^{N}{\left( {x_{k} - \mu_{x}} \right)\left( {y_{k} - \mu_{y}} \right)}}{N - 1}} & (27)\end{matrix}$where C_(i,j) is the covariance of variable i with variable j, x_(k) andy_(k) are the ith and jth feature values of the kth object, μ_(x) is themean of all N values of x_(k), μ_(y) is the mean of all N values ofy_(k). The covariance matrix is 3×3 and PCA yields three eigenvectors:the eigenvectors are treated as a 3×3 matrix, which is used to multiplythe transpose of the N×3 image matrix to produce a product matrix. Theproduct matrix has an N×1 first column which is the first principalcomponent, which may be considered as the most important component. Itis the component with the maximum eigenvalue, and it provides agreyscale sub-image (one pixel value for each of N pixels) with amaximum range of information compared to equivalents associated withother components. PCA is carried out for each of the overlapping windowsdefined above and each provides a respective first principal componentand greyscale sub-image of size 128×128 pixels.

At 122, each sub-image resulting from 120 is converted to acorresponding binary sub-image by applying a thresholding methodreferred to as “Otsu”. Otsu is a standard thresholding techniquepublished by Otsu N., ‘A thresholding selection method from greylevel-histograms’, IEEE Trans Systems, Man & Cybernetics, vol. 9, 1979,pp 62-66. The Otsu threshold selection method aims to minimise for twoclasses a ratio of between-class variance to within-class variance: i.e.the higher the variance between classes the better the separation. Inthe present example the two classes are a below-threshold class (pixelvalue 0) and an above-threshold class (pixel value 1), so by applyingOtsu thresholding the greyscale sub-image is converted into a binarysub-image containing a set of blobs.

At 124 all blobs (objects) that touch or intersect sub-image boundariesare removed. Thus, if at any pixel a blob meets a border it is removedby setting its pixels to a background pixel value. This is because suchboundaries give blobs meeting them artificial straight edges which cangive misleading results later. Because of sub-image overlap, a blobwhich appears partly in one image may appear wholly in anothersub-image. This step 124 only arises from the use of sub-images.

At 126 the outputs from 124 are inverted and connected componentlabelling (CCL, as described earlier) is applied in order to enable anyholes in the blobs to be removed. This is not essential, but it providesspatial filtering which improves results somewhat. Because of theinversion, areas of pixel value 1 labelled by CCL labelling will bebackground pixels and holes within blobs. Holes, i.e. all labelled areasother than background pixels, are removed (filled in) by setting pixelsof holes in each blob to the value of the blob's other pixels.

At 128 the outputs from 126 are inverted once more and CCL is applied:after this inversion, the labelled areas are the blobs within thesub-image filled at 126. CCL also yields blob centre positions usedlater. Any blobs smaller than a minimum area threshold of 400 pixels arerejected, i.e. set to an image background value, in accordance with thetable immediately below. This is another desirable but not essentialspatial filtering step. Blob size < min area (400 pixels) Reject blobOtherwise Accept for further processing into set of labelled blobs

At 130 multiple sub-images output at 128 are reassembled into a newbinary image which is the same size as the original before decompositionat 120; the new binary image has undergone filtering and now containsonly blobs that are of interest for subsequent mitosis processing. Imagepreprocessing terminates with step 130: the set of blobs remaining hasbeen cleared both of unwanted small blobs and of holes within blobs.Image preprocessing using steps 120-130 is advantageous because it doesnot significantly affect shapes of blob perimeters, which is importantfor mitosis analysis. Results from whichever of the mitosis cueingprocesses 16 and 18 has been used are introduced in this step 130: i.e.a blob is accepted or rejected according to whether or not there is ablob in the mitosis cued image in substantially the same position (thiscould be implemented lo using a logical AND operation). It is notessential to use a process 16 or 18 here but it is useful to reduceprocessing burden.

At 132 Principal Component Analysis (PCA, as previously described) isapplied to the entire input (RGB) image. As before, ae red, green orblue image plane or an average thereof could be used, but PCA ispreferred. PCA yields a first component which is a greyscale image witha better information range than those of other components. A featureextraction procedure 134 is applied to a local window of 51×51 pixelscentered on the centres of each of the blobs identified inpre-processing at 130 and appearing in the greyscale image. Theprocedure 134 determines an average cross-section (profile) for thatrespective region of the greyscale image that corresponds to each blob:for the purposes of this calculation each grayscale value used isnormalised to lie in the range 0 to 1 by dividing it by 255. Arespective series of profiles of each blob is taken using a linefifty-one pixels long extending across the respective greyscale imageregion corresponding to that blob and centered on the blob centre: thisgives fifty-one pixel values or histogram points per profile, andprofiles are taken at nine different angular orientations at 20 degreeintervals: the respective mean of the nine profiles of each blob is thencalculated.

A respective histogram of each mean profile is then obtained andquantised to just five intervals or bins 1 to 5 as follows, (1)0≦profile<0.2, (2) 0.2≦profile<0.4, (3) 0.4≦profile<0.6, (4)0.6≦profile<0.8, and (5) 0.8≦profile≦1.0. The bins have centres at 0.1,0.3, 0.5, 0.7 and 0.9 respectively. Each bin contains the number ofpixels in the mean profile having its respective intensity value:because this number is averaged over the nine profiles it need not be aninteger. Bin 1 corresponds to a darkest group of image intensity values,i.e. low greyscale values of the kind one would associate with images ofmitotic figures; mitotic figures are normally dark using conventionalstaining techniques so relatively darker degrees of grey level are ofmore interest bins 2 to 4 correspond to progressively brighter valuesand bin 5 to the brightest of the five values. These are relative thoughbecause the profiles all come from relatively dark image regions. Anapproximate mean profile is represented by five values each of which isan intensity value in a respective bin averaged over nine measuredprofiles. Each set of five values characterises a blob now treated asindicating an actual cell. The minimum value of each mean profile isrecorded as the variable ‘minprofile’, this being the contents (numberof pixels averaged over nine profiles) of the bin having the smallestcontents of all five bins.

At 136 the contents of bins 1, 2 and 4 are used to determine if acurrent cell that corresponds to a current blob is mitotic or not.Specifically, the following criteria are applied for each blob, where“bin (n)” means the contents of the nth bin and n=1, 2 or 4:

If bin (1)>7.6, then:

-   -   -   current cell is mitotic.

Otherwise, if bin (1)≦7.6, then

-   -   if bin (2)>25.5 and bin (4)>0 and minprofile<0.15, then:        -   current cell is mitotic    -   otherwise, i.e. if bin (2)≦25.5 and/or bin (4)≦0 and/or        minprofile≧0.15) then        -   current cell is not mitotic

The first criterion—bin (1)>7.6—may if desired be the only one used todetermine whether a cell is mitotic. The “otherwise” criterion—bin(1)≦7.6—is optional.

Each alternative mitotic feature detection technique 20 and 24 producesmeasurements derived from ten images. Each mitotic feature detectiontechnique is applied to ten images or tiles as has been said: themitotic figures are counted for each image and the counts are addedtogether to provide a total for the ten images. The mitotic figure countfor a technique is low, medium or high with points 1, 2 or 3 accordingto whether it is 0 to 5, 6 to 10 or 11 or more respectively as shown inthe table below. Measurement: Number of Mitotic figures in Ten ImagesMeaning Points 0-5 Low 1  6-10 Moderate 2 ≧11 High 3

The measurement of mitosis may be combined with others obtained forpleomorphism and tubules by different methods to derive an overallgrading referred to in medicine as a “Bloom and Richardson grading”: itis used by clinicians as a measure of cancer status.

The examples given in the foregoing description for calculatingintermediate quantities and results can clearly be evaluated by anappropriate computer program recorded on a carrier medium and running ona conventional computer system. Examples of program steps haven beengiven. Such a program is straightforward for a skilled programmer toimplement without requiring invention, because the procedures are wellknown. Such a program and system will therefore not be describedfurther.

1. A method of measurement of mitotic activity from histopathologicalspecimen image data, the method comprising the steps of: a) identifyingpixels in the image data having luminances associated with mitoticfigures; b) selecting from among the identified pixels a reference pixelwhich is sufficiently close in position and luminance to anotheridentified pixel to provide a reference colour; c) locating pixels inthe image data with luminances sufficiently close to that of thereference colour to indicate potentially mitotic figures; d)incrementing image regions corresponding to potentially mitotic figuresfrom the located pixels by adding pixels thereto, potential incrementsto image regions being implemented or rejected by according to whetheror not their luminances are sufficiently close to respective imageregion luminances and sufficiently far from an image data backgroundluminance; e) selecting grown image regions on the basis of thresholdsfor image region area, compactness and width/height ratio; and f)counting selected grown image regions as actually indicating mitoticfigures on the basis of a threshold for number of such regions.
 2. Amethod according to claim 1 wherein the step of selecting grown imageregions also involves thresholds for ratio of image region luminance tobackground luminance and area difference between areas derived bygrowing each image region with multiple thresholds.
 3. A methodaccording to claim 2 wherein the thresholds for image region area,compactness, width/height ratio, luminance and area difference are: 355pixels<area<1700 pixels, 0.17<compactness<0.77, width/height ratio<2.7,luminance percentage<44%, area difference<23area/100.
 4. A methodaccording to claim 1 wherein the step of counting selected grown imageregions as actually indicating mitotic figures also involves thresholdsfor region area and luminance.
 5. A method according to claim 1 whereinsuccessive potential increments to image regions are individual pixelseach of which is an immediate row or column neighbour of an existingimage region pixel.
 6. A method according to claim 1 wherein step b) isimplemented with a reference pixel having a luminance differing by lessthan about 8% compared to another identified pixel distant from it bynot more than two percent of a smaller of two image dimensions.
 7. Amethod according to claim 1 wherein step a) includes white balancing andmedian filtering the image data prior to identifying pixels havingluminances corresponding to mitotic figures.
 8. A method according toclaim 1 wherein in step c) pixels are cued for acceptance or rejectionas regards indicating mitotic figures by: a) thresholding colour imagedata to remove pixels lacking intensities associated with mitotic figureimagery, b) removal pixels not present in all colours, and c)thresholding image region areas to remove those too small and too largeto be potential mitotic figures.
 9. A method according to claim 1wherein in step c) pixels are cued for acceptance or rejection asregards indicating mitotic figures by: a) segmenting to identify pixelswith intensities associated with mitotic figure imagery, b) thresholdingimage region areas to remove those too small and too large to bepotential mitotic figures, c) cluster analysis to determine whether ornot a pixel's image region is in a sufficiently large cluster, and d)necrotic and hairy edge filtering.
 10. A method of measuring mitoticactivity from histopathological specimen image data, the method havingthe steps of: a) measuring an intensity profile of an image regioncorresponding to a potentially mitotic figure, and b) counting the imageregion as indicating a mitotic figure if its profile has a value greaterthan a prearranged threshold at a position in the profile havingintensity associated with mitotic figure imagery.
 11. A method accordingto claim 10 including counting the image region as indicating a mitoticfigure if its profile has a first value not greater than the prearrangedthreshold at a position in the profile having intensity associated withmitotic figure imagery, a second value greater than a prearranged secondthreshold, a third value greater than a prearranged third threshold, anda minimum value less than a prearranged fourth threshold.
 12. A methodaccording to claim 11 wherein the first value is at one end of theprofile, the first and second values adjoin one another in the profileand the third value does not adjoin the second value.
 13. A methodaccording to claim 11 wherein the image data comprise a first PrincipalComponent obtained by Principal Component Analysis (PCA) of colouredimage data.
 14. A method according to claim 11 wherein step a) includespreprocessing image data by: a) decomposing the image data intooverlapping sub-images, b) applying PCA to the sub-images to derive afirst Principal Component image, c) thresholding the first PrincipalComponent image to produce a binary image of blobs and background d)rejecting blobs adjacent to or intersecting sub-image boundaries, e)filling holes in blobs, f) rejecting blobs too small to correspond topotential mitotic figures, and g) reassembling the sub-images into asingle image for image region profile measurement as aforesaid in stepa).
 15. A method according to claim 14 wherein after step g) pixels arecued for acceptance or rejection as regards indicating mitotic figuresby: a) thresholding colour image data to remove pixels lackingintensities associated with mitotic figure imagery, b) removal pixelsnot present in all colours, and c) thresholding image region areas toremove those too small and too large to be potential mitotic figures.16. A method according to claim 14 wherein after step g) pixels are cuedfor acceptance or rejection as regards indicating mitotic figures by: a)segmenting to identify pixels with intensities associated with mitoticfigure imagery, b) thresholding image region areas to remove those toosmall and too large to be potential mitotic figures, c) cluster analysisto determine whether or not a pixel's image region is in a sufficientlylarge cluster, and d) necrotic and hairy edge filtering.
 17. Computerapparatus for measuring mitotic activity from histopathological specimenimage data, the apparatus being programmed to execute the steps of: a)identifying pixels in the image data having luminances associated withmitotic figures; b) selecting from among the identified pixels areference pixel which is sufficiently close in position and luminance toanother identified pixel to provide a reference colour; c) locatingpixels in the image data with luminances sufficiently close to that ofthe reference colour to indicate potentially mitotic figures; d)incrementing image regions corresponding to potentially mitotic figuresfrom the located pixels by adding pixels thereto, potential incrementsto image regions being implemented or rejected by according to whetheror not their luminances are sufficiently close to respective imageregion luminances and sufficiently far from an image data backgroundluminance; e) selecting grown image regions on the basis of thresholdsfor image region area, compactness and width/height ratio; and f)counting selected grown image regions as actually indicating mitoticfigures on the basis of a threshold for number of such regions. 18.Apparatus according to claim 17 programmed to execute the step ofselecting grown image regions by also using thresholds for ratio ofimage region luminance to background luminance and area differencebetween areas derived by growing each image region with multiplethresholds.
 19. Apparatus according to claim 18 wherein the thresholdsfor image region area, compactness, width/height ratio, luminance andarea difference are: 355 pixels<area<1700 pixels, 0.17<compactness<0.77,width/height ratio<2.7, luminance percentage<44%, areadifference<23area/100.
 20. Apparatus according to claim 17 programmed toexecute the step of counting selected grown image regions as actuallyindicating mitotic figures by also using thresholds for region area andluminance.
 21. Apparatus according to claim 17 wherein successivepotential increments to image regions are individual pixels each ofwhich is an immediate row or column neighbour of an existing imageregion pixel.
 22. Apparatus according to claim 17 programmed to executestep b) with a reference pixel having a luminance differing by less than8% compared to another identified pixel distant from it by not more thantwo percent of a smaller of two image dimensions.
 23. Computer apparatusfor measuring mitotic activity from histopathological specimen imagedata, the apparatus being programmed to execute the steps of: a)measuring an intensity profile of an image region corresponding to apotentially mitotic figure, and b) counting the image region asindicating a mitotic figure if its profile has a value greater than aprearranged threshold at a position in the profile having intensityassociated with mitotic figure imagery.
 24. Apparatus according to claim23 programmed to count an image region as indicating a mitotic figure ifits profile has a first value not greater than the prearranged thresholdat a position in the profile having intensity associated with mitoticfigure imagery, a second value greater than a prearranged secondthreshold, a third value greater than a prearranged third threshold, anda minimum value less than a prearranged fourth threshold.
 25. Apparatusaccording to claim 24 wherein the first value is at one end of theprofile, the first and second values adjoin one another in the profileand the third value does not adjoin the second value.
 26. Apparatusaccording to claim 24 wherein the image data comprise a first PrincipalComponent obtained by Principal Component Analysis (PCA) of colouredimage data.
 27. Computer program code for use in measuring mitoticactivity from histopathological specimen image data, the computerprogram code containing instructions to control a computer to implementthe steps of: a) identifying pixels in the image data having luminancesassociated with mitotic figures; b) selecting from among the identifiedpixels a reference pixel which is sufficiently close in position andluminance to another identified pixel to provide a reference colour; c)locating pixels in the image data with luminances sufficiently close tothat of the reference colour to indicate potentially mitotic figures; d)incrementing image regions corresponding to potentially mitotic figuresfrom the located pixels by adding pixels thereto, potential incrementsto image regions being implemented or rejected by according to whetheror not their luminances are sufficiently close to respective imageregion luminances and sufficiently far from an image data backgroundluminance; e) selecting grown image regions on the basis of thresholdsfor image region area, compactness and width/height ratio; and f)counting selected grown image regions as actually indicating mitoticfigures on the basis of a threshold for number of such regions. 28.Computer program code according to claim 27 incorporating instructionsto provide for implementing the step of selecting grown image regions byalso using thresholds for ratio of image region luminance to backgroundluminance and area difference between areas derived by growing eachimage region with multiple thresholds.
 29. Computer program codeaccording to claim 28 wherein the thresholds for image region area,compactness, width/height ratio, luminance and area difference are: 355pixels<area<1700 pixels, 0.17<compactness<0.77, width/height ratio<2.7,luminance percentage<44%, area difference<23area/100.
 30. Computerprogram code according to claim 27 incorporating instructions to providefor implementing the step of counting selected grown image regions asactually indicating mitotic figures using also thresholds for regionarea and luminance.
 31. Computer program code for use in measuringmitotic activity from histopathological specimen image data,incorporating instructions to provide for implementing the steps of: a)measuring an intensity profile of an image region corresponding to apotentially mitotic figure, and b) counting the image region asindicating a mitotic figure if its profile has a value greater than aprearranged threshold at a position in the profile having intensityassociated with mitotic figure imagery.
 32. A computer program accordingto claim 31 incorporating instructions to provide for counting the imageregion as indicating a mitotic figure if its profile has a first valuenot greater than the prearranged threshold at a position in the profilehaving intensity associated with mitotic figure imagery, a second valuegreater than a prearranged second threshold, a third value greater thana prearranged third threshold, and a minimum value less than aprearranged fourth threshold.
 33. A computer program according to claim32 wherein the first value is at one end of the profile, the first andsecond values adjoin one another in the profile and the third value doesnot adjoin the second value.
 34. A computer program according to claim32 incorporating instructions to provide for step a) to includepreprocessing image data by: a) decomposing the image data intooverlapping sub-images, b) applying PCA to the sub-images to derive afirst Principal Component image, c) thresholding the first PrincipalComponent image to produce a binary image of blobs and background d)rejecting blobs adjacent to or intersecting sub-image boundaries, e)filling holes in blobs, f) rejecting blobs too small to correspond topotential mitotic figures, and g) reassembling the sub-images into asingle image for image region profile measurement as aforesaid in stepa).